#include <stdio.h>
#include <string.h>
#include <math.h>
#include <ctype.h>
#include <stdlib.h>
#define ll long long
#define maxn 100+5
#pragma GCC optimize(2)
const double eps=1e-15;

int powerMethord1(double a[][maxn],double y[][maxn],double beta[],int n);
int powerMethord2(double a[][maxn],double y[][maxn],double beta[],int n);
double dotProduct(double a[],double b[],int n);
void martixVectorMal(double a[][maxn],double x[],double b[],int n);
double Abs(double x);
void swap(int *x,int *y);
double infinityNorm(double x[],int n);
double sgn(double x);


int main(){
    double a[maxn][maxn],y[2][maxn],beta[2];
    int n;
    
    printf("输入矩阵阶数n:");
    scanf("%d",&n);
    printf("输入%d阶矩阵:\n",n);
    for(int i=1;i<=n;++i)
        for(int j=1;j<=n;++j)
            scanf("%lf",&a[i][j]);
    printf("选择迭代格式：\n1:2-范数迭代\n2:无穷范数迭代\n");
    int op;
    scanf("%d",&op);
    int now=op==1?powerMethord1(a,y,beta,n):powerMethord2(a,y,beta,n);
    if(now>=0){
        printf("绝对值最大的特征值为%10.6lf\n",beta[now]);
        printf("对应的特征向量为:[");
        for(int i=1;i<=n;++i)
            printf("%10.6lf",y[now==1?0:1][i]);
        printf("]");
    } else 
        printf("迭代失败,请检查收敛性。\n");
    return 0;
}

inline int powerMethord1(double a[][maxn],double y[][maxn],double beta[],int n){
    double u[2][maxn],eta[2];
    int count=0,now=1,pre=0;
    memset(u,0,sizeof(u));
    u[pre][n]=1.0;
    while(count<=10000000){
        eta[pre]=sqrt(dotProduct(u[pre],u[pre],n));
        for(int i=1;i<=n;++i)
            y[pre][i]=u[pre][i]/eta[pre];
        martixVectorMal(a,y[pre],u[now],n);
        beta[now]=dotProduct(y[pre],u[now],n);
        if(Abs(beta[now]-beta[pre])/Abs(beta[now]) < eps){
            printf("共迭代%d次\n",count+1);
            return now;
        } else {
            swap(&pre,&now);
        }
        ++count;
    }
    return -1;
}

inline int powerMethord2(double a[][maxn],double y[][maxn],double beta[],int n){
    double u[2][maxn];
    int count=0,now=1,pre=0;
    memset(u,0,sizeof(u));
    u[pre][n]=1.0;
    while(count<=10000000){
        double h[2],maxAbs=-1.0;
        h[pre]=infinityNorm(u[pre],n);
        
        for(int i=1;i<=n;++i)
            y[pre][i]=u[pre][i]/Abs(h[pre]);
        martixVectorMal(a,y[pre],u[now],n);

        h[now]=infinityNorm(u[now],n);
        beta[now]=sgn(h[pre])*h[now];

        if(Abs(beta[now]-beta[pre])/Abs(beta[now]) < eps){
            printf("共迭代%d次\n",count+1);
            return now;
        } else {
            swap(&pre,&now);
        }
        ++count;
    }
    return -1;
}

double dotProduct(double a[],double b[],int n){
    double sum=0.0;
    for(int i=1;i<=n;++i)
        sum+=a[i]*b[i];
    return sum;
}

inline void martixVectorMal(double a[][maxn],double x[],double b[],int n){
    for(int i=1;i<=n;++i){
        double sum=0.0;
        for(int j=1;j<=n;++j)
            sum+=a[i][j]*x[j];
        b[i]=sum;
    }
}

inline double infinityNorm(double x[],int n){
    double ans,maxAbs=-1.0;
    for(int i=1;i<=n;++i){
        if(Abs(x[i])>=maxAbs){
            maxAbs=Abs(x[i]);
            ans=x[i];
        }
    }
    return ans;
}

inline double Abs(double x){
    return x<0.0?-x:x;
}

inline void swap(int *x,int *y){
    int temp=*x;
    *x=*y;
    *y=temp;
}

inline double sgn(double x){
    return x>0.0?1.0:(x<0.0?-1:0.0);
}